Multiple co-frequency sources DOA estimation for coprime vector sensor arrays

For the problem of direction-of-arrival (DOA) estimation using a coprime array, there are high spatial spectrum outputs of false alarms caused by the overlap of main and grating lobes from subarrays. In this paper, a DOA estimation method of more than two co-frequency sources for a coprime vector hydrophone array is proposed. The method is based on vector cross terms (VCTs), making full use of the directivity of channel combinations for vector hydrophones. Based on VCTs, the characteristic data point identification method is conducted and ensures that the bearing data with the characteristic can be preserved. For further interference rejection, the paper designs Queue Selection (QS) method based on inverse beamforming. The influence of grating lobes can be weakened with the QS, further improving the accuracy of direction extraction. The algorithm in this work does not require decoherence processing, and the simulation work shows that it achieves stable DOA estimation with a low signal-to-noise ratio (SNR).


Introduction
Coprime linear arrays (CLAs) have become popular nowadays for providing system array configurations beyond Nyquist sampling and reducing mutual coupling between array elements [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15][16][17]. The classical spatial spectrum estimation algorithms for CLAs are based on Product [18] and Min [19] theories. The Product theory calculates the spatial power spectral density (PSD) of the observed signal through multiplication. The multiplication is conducted by a beamformed subarray scanned response with the complex conjugate of another. The Min processor estimates the spatial PSD by selecting the lowest value from the CLA subarray periodograms for each bearing. The Product method can resolve the spatial aliasing ambiguity arising from CLA subarrays, and the Min processor maintains identical array resolution while attaining the reduced height of peak sidelobe and total sidelobe area with the same extended CLA geometry.
Direction of arrival (DOA) estimation is a vital area of array signal processing in acoustics, radar, and wireless communications [20]. For the issue of DOA estimation of CLAs, the DOA estimation method obtained by combining the multiple signal classification (MUSIC) results a1111111111 a1111111111 a1111111111 a1111111111 a1111111111 of the linear arrays (DECOM) utilizes a combination of the MUSIC results from two decomposed uniform linear arrays to obtain DOA estimates [21]. An underwater vector sensor coprime line array (vCLA) with an extended configuration and sound pressure and vibration velocity joint processing (PVJP) was proposed for improving beamforming and DOA estimation performance [22]. The DECOM method and the extended vCLA are based on the physical structure of the array itself, as opposed to some methods, based on research conducted in the virtual array element domain. For example, a new DOA estimation algorithm utilizing Toeplitz matrix reconstruction and virtual array interpolation methods was introduced [23]. A discrete inverse Fourier transform (IDFT) method employed virtual signals of second order with angular-spatial frequencies to estimate the DOA information [24]. These methods focus on the advantage of high degrees of freedom (DOFs) using the virtual array elements and are primarily used in the radar field. In [25], a technique was presented that utilized a coprime electromagnetic vector-sensor (EMVS) array. The method relied on solving a nuclear norm minimization (NNM) problem to construct a larger covariance matrix that was then used for DOA estimation. Then the problem of DOA estimation in two dimensions was discussed for a bistatic multiple-input multiple-output (MIMO) radar system with coprime EMVS arrays in [26]. There is also some research about DOA estimation in different scenarios for virtual array domain processing. A multi-source DOA tracking method was presented [27]. The number of sources was estimated using the minimum description length (MDL), and then the MUSIC pseudo-spectrum was utilized. For impulsive noise and non-circular signals scenarios, the DOA estimation algorithm for a coprime array was shown [28], which proposed an augmented phased fractional low-order moment (A-PFLOM) to suppress noise and applied reduced-dimension MUSIC (RD-MUSIC) subspace techniques for DOA information. For mixed noise scenarios, a DOA estimation method based on phased fractional low-order moment was advanced [29]. For the scenario of impulsive noise, fractional low-order moments were applied [30].
Nevertheless, CLAs suffer from significant limitations in practical applications. High grating and peak side lobes degrade the performance of some spatial spectrum estimation algorithms and even cause many false alarms in target detection. Many scholars have conducted research to improve this problem. An extended coprime array method was proposed to improve the peak side lobe of the CLAs. The method derived a sufficient condition on the extension factor connected with the location and formation of peak side lobes [6]. Similarly, the approaches based on extended CLAs were developed as semi-coprime arrays [10], threelevel and multi-level prime arrays [7,8,11,12]. Additionally, designing an array factor by nulling the grating lobes with known directions can also be a solution for this issue [31]. The above problems and solution ideas are based on incoherent signals, while coherent signals are less considered. For coherent signals, especially co-frequency signals, there are more severe problems. The spatial spectrum output of the co-frequency signals of the CLAs will produce grating lobes equivalent to the spatial spectrum of the targets' directions, which undoubtedly severely impacts the accurate detection of the targets [14]. High-resolution DOA estimation algorithms like MUSIC [32] and estimating signal parameters via rotational invariance techniques (ESPRIT) [33] will fail because of the source covariance matrix's rank deficiency when the incident signals are coherent. For coherent target detection, several methods have been proposed, including spatial smoothing (SS) [34] and forward/backward SS (FBSS) [35]. However, the SS technique sacrifices array aperture to attain DOA estimation of coherent signals, while the FBSS approach yields superior estimation precision but the incomplete implementation of signal decoherence. In addition, uniform linear arrays (ULAs) are often given more attention than sparse linear arrays in these methods.
For underwater acoustic signal detection and estimation, the grating lobe effect of the array can be reflected in sparse arrays. For instance, a vector hydrophone linear array performs similarly to a sound pressure array for detecting and estimating underwater acoustic signals, but with a smaller array aperture [36]. Furthermore, the sound pressure and vibration velocity combination exhibits robust anti-isotropic noise capabilities [37][38][39]. In practical applications, vector hydrophone linear arrays allow for a certain degree of sparsity, but the accompanying grating lobes issues are still expected to be addressed.
In this paper, we develop a DOA estimation algorithm for more than two co-frequency signals using a coprime vector sensor array. The conventional beamformer for the whole coprime array is used as the pre-processor, and the vector cross terms (VCTs) are constructed on two subarrays. The algorithm uses the channel combination of vector hydrophones, and the characteristic data point identification method based on VCTs is designed. In order to improve the accuracy, a queue selection of signals based on inverse beamforming is used to perform filtering of suspicious targets. Finally, the orientation information of real targets can be obtained. The proposed method differs from existing techniques in handling co-frequency signals without spatial smoothing, thereby reducing redundancy. The simulation results demonstrate the effectiveness of the proposed algorithm.
Notations: The matrices (vectors) are represented by upper-case (lowercase) bold characters. (�)*, (�) H and (�) T denote the complex conjugate, conjugate transpose and transpose, respectively. I stands for the unit matrix. � represents the Kronecker product.

Coprime vector sensor array signal model
An underwater acoustic vector sensor coprime linear array (vCLA) is made up of two sparse uniform subarrays of vector sensors. The first subarray has M physical sensors, while the second has N physical sensors. It is important to note that M and N are coprime numbers. In Fig 1, it can be shown that the first subarray, comprising M sensors, has an intersensor distance of Nd, whereas the second subarray, consisting of N sensors, has an intersensor spacing of Md.
Here, d = λ/2 refers to the unit array element spacing, where λ represents the wavelength of the signal impinging on the array. Then the representation of the array configuration is as follows.
Supposing that the desired signal comes from direction θ 0 , then we can model the received signal as follows. where T refers to the vector of the signal waveform and NðtÞ ¼ ½n T 1 ðtÞ; n T 2 ðtÞ; � � � n T MþNÀ 1 ðtÞ� T � CN ð0; s 2 n IÞ represents the Gaussian noise component with statistical independence and the power of s 2 n . Here n i (t) = [n p (t), n vx (t), n vy (t)] T , i = 1, 2, � � � M + N − 1 is the noise vector, which is composed of the pressure component and the horizontal velocity components, denoted as n p (t), n vx (t) and n vy (t) respectively. a(θ) is the steering vector connected with DOA θ 0 : The acoustic vector sensor is a device that can simultaneously measure sound pressure and particle velocity at a single point in the sound field, reducing the need for multiple devices and minimizing measurement redundancy. As shown in Fig 2, the relationship between different channels of a vector sensor after rotation can be obtained from the Givens transform as follows: So the 3 × 1 steering vector can be formulated as: The conventional beamformer (CBF) is given by [40]: where w CLA is the weight of beamformer and d CLA 2 S. x vCLA denotes the signal vector of the vCLA. In practical applications, the maximum beam output in the direction of the source can be achieved when there is a single source. Nevertheless, When sources have the same frequency, the significant cross term arises in Eq (6). Consequently, there will be high beam outputs in directions other than the target directions, which could cause false alarms or inaccurate estimations of the DOAs. The cross term can be attributed to the relative positioning of the main lobes and grating lobes in the beam output of two subarrays in a coprime vector sensor array [41]. There is the following relationship along the axis of sin(θ).
where i 1 = 0, 1, 2, . . ., i 2 = 0, 1, 2, . . ., and θ N and θ M are two angles causing false alarms. When the main lobe of one subarray aligns with the grating lobe of the other subarray, or when the grating lobe of one subarray aligns with the grating lobe of the other subarray, the output beams have similar amplitudes to the true sources, leading to ambiguity.

DOA estimation for multiple coherent sources with same frequency
The VCTs construction for coprime vector hydrophone array The proper combination of sound pressure and vibration velocities can form various directivity of a vector hydrophone. A vector hydrophone can effectively reduce noise by combining

PLOS ONE
sound pressure and vibration velocity in the processing of the acoustic vector signal. Therefore, to obtain the desired output, the channels of acoustic vector hydrophones are transformed by rotation and combination methods. The Eq (10) is obtained Here, φ denotes the angle of electron rotation and S(t) represents the sound pressure signal. This paper utilizes a combination of sound pressure and vibration velocity: where It can be seen that B s (θ) = 0 when θ = φ. The method used for noise reduction by combining sound pressure and vibration velocity processing can be classified as a generalized spatial domain filtering technique. Therefore, selecting a suitable rotation angle φ and rotating the acoustic vector hydrophone data can reduce noise and effectively detect weak targets. For a single uniform sparse vector hydrophone array, traditional beamforming based on Eq (13) can generate a spatial spectrum output related to the target direction (take the target from 70˚as an example). The output is shown in Fig 3, and "Normal Channels" in the legend represents the sound pressure and velocity channels of vector hydrophones are treated as the same type channels, i.e., the data of each channel is not combined. Without considering noise, the spatial spectrum output based on a vector hydrophone combined channels exhibits a concave point at the source direction, which is shown in Fig 3(a). Unfortunately, the reliability of using concave points to determine target directions may be affected by the presence of noise. What's worse, for a sparse array, the unreliability can be more severe as the array spacing increases. In this paper, we proposed a DOA estimation method for a coprime vector hydrophone array, which uses vector cross terms (VCTs) inspired by the Product theorem [18,41]. This method resolves spatial frequency ambiguities by multiplying the response of one subarray with the complex conjugate of the other [1]. Suppose that the vector hydrophones in a subarray consisting of M sensors receive the data of the acoustic pressure, x-axial acoustic particle velocity, and y-axial acoustic particle velocity, which can be denoted as P M (t), V xM (t), and V yM (t), respectively. On the other hand, let P N (t), V xN (t), and V yN (t) denote the data of acoustic pressure, x-axial acoustic particle velocity, and y-axial acoustic particle velocity that the vector hydrophones from a subarray with N sensors have received. Then we can create the vector cross terms for a coprime vector hydrophone array by using the following construction: where

Data cross validation based on VCTs
The concave point discriminant function we introduce is based on the output of the array spatial spectrum and is designed to identify concave points.
where b NM (θ) and b MN (θ) are the beam outputs calculated by using the covariance of combined data generated by Eq (14). The reliability of concave point judgment can be improved by establishing a relationship between two sparse subarrays of the coprime vector hydrophone array instead of using a single sparse array. The inverse Advance and Retreat method [42] is utilized to detect concave points in the direction of suspected targets. Let F be the search step and ϕ s be the suspected target's orientation. The discriminating process can be expressed as � IFðpð� s Þ À pð� s þ FÞ < 0Þ; � IFðpð� s þ FÞ À pð� s þ 2FÞ < 0Þ; D p;À F ¼ IFðpð� s À FÞ À pð� s À 2FÞ < 0Þ � IFðpð� s À FÞ À pð� s Þ < 0Þ: where 'IF()' represents an if conditional operation. The algorithm diagram for the discriminating process is shown as Fig 4, where � 0 r stands for the orientation result of discriminating. The discriminating process of each orientation point is based on the joint judgment of its surrounding expansion points, and the minimum interval is related to the algorithm resolution.

Queue Selection (QS) method of multi-frequency coherent targets based on signal energy
For the case of two targets with the same frequency, the method mentioned above can directly obtain the directions of real targets. Nevertheless, when the number of targets exceeds two, the interference between the grating lobes caused by the same frequency signal intensifies, resulting in the inability to obtain the true DOAs directly. In order to reduce this kind of impact, this article uses inverse beamforming cancellation technology to select the estimated direction of arrival one by one in the queue, thereby reducing false alarms caused by interference between co-frequency grating lobes [43]. Let Θ be a set of predicted directions.
where S 1 denotes the number of suspected targets obtained based on Eq (18) in the DOA preestimation step. For example, assume that an interference target comes in the direction of θ 1 . For the subarray with M sensors of the coprime vector sensor array, the array element response vector is expressed as:

PLOS ONE
interference bearing is located can be expressed as: where x(t) is the actual signal. The sector is centered on θ 1 and contains 2O + 1 beams. The array weighting matrix corresponding to the b-th beam is formulated as: where φ bM = φ M . By inverse beamforming, we obtain: When the sector is approximately a small angle sector, e {−j[(i−1)(φ−φbM)]} ffi 1, then I ωM = s m .
Then we can derive the information for the velocity of the vector array: Similarly, the operation formulas for the subarray with N sensors are as follows: a N ðy 1 Þ ¼ ½1; e À jφ N ; e À j2φ N ; . . . ; e À jðNÀ 1Þφ N � ð28Þ where φ bN = φ N . The sound pressure and velocity channels are used as general channels to obtain beam data in the interference direction θ 1 . Let the signal reception matrix after interference rejection be Z, then the received data of the incident signal after cancellation is: In the process of interference rejection, the grating lobes brought by the interference are also weakened, which plays a role in the discrimination of whether the target is real. Queue interference suppression is performed in the pre-estimated orientation target obtained based on Eq (10) to Eq (18), and then judge whether the pre-estimated directions still exist in the result, so finally, true targets orientations can be obtained. The pseudo-code of this part is presented in Algorithm 1.

Algorithm 1 Pseudo code of the QS method.
Initailize: 1: Input data: Array signal S, V x and V y , Interference direction θ i 2: Initailize parameters: Signal integral length T s , Angle search range Θ s Main Loop: 3: while Length(S t ) do = T s 4: Array element response vector a(θ i ) with Eq (20) and Eq (28) 5: Inverse beamforming weighting matrix W with Eq (24) and Eq (33) 6: Array beam output Z with Eq (25) to Eq (27) and Eq (34) to Eq (36) 7: Update Θ c with Eq (14) to Eq (18) 8: Array data after suppression S*, V * x and V * y 9: for i = 1 : Θ s do 10: Beamforming B i with Eq (6) 11: end for 12: end while 13: Output B Y s For example, as shown in Fig 5, there are two sources and one interference. The directions of two sources are 50˚and 80˚respectively and the interference comes from 120˚. For the sparse vector sensor array, multiple targets and interferences and their grating lobes are in the array spatial spectrum output before interference suppression. After the interference of 120˚is suppressed, the grating lobe effect produced by it also disappears, and then the grating lobe suppressing is achieved.

Major steps and practical application
(1) First, the number of suspected targets and their respective orientations are determined in advance based on Eq (6).
(2) Second, the VCTs and discriminant function to identify concave points are created using Eqs (14) and (17).
(3) The direction of the source can be determined in both scenarios where a single target is present or when two targets are detected with a defined detection threshold. As there is no chance of false alarms in the array's output in both scenarios, only the actual output is displayed. The target DOA estimation preprocessing can be achieved by Eqs (9) and (18). For more directions of multiple co-frequency signals, the queue selection processing based on Eqs (19) to (39) can be conducted for interference rejection.
(4) Then interference suppression and target orientation verification are performed sequentially to obtain the final target orientation results, which are executed by Eqs (17) to (39). The pseudo code of the overall algorithm is shown in Algorithm 2.
Algorithm 2 Pseudo code of the major steps for the overall algorithm.

Initailize:
Input data: Array signal S t

Simulation and results
Several scenarios were simulated in this part. In the simulation work, a coprime array with M = 5 and N = 6 sensors for two subarrays, respectively. The number of elements for the whole array is M + N − 1 = 10. Comparisons between the proposed method and MUSIC based on the FBSS algorithm are included. Firstly, when there are three sources from the directions of 50˚, 70˚, and 110˚respectively, with SNR = 0 dB and f=500 Hz. As shown in Fig 6(a), when the three co-frequency targets' directions do not hold the grating lobes relationship, the array can directly output the accurate orientations without false alarm interferences. The results from both algorithms are consistent with the directions of the actual source. In Fig 6(b), it can be seen that if two of the three targets have a grating lobe relationship and the third target is not in the relationship of the grating lobes with the first two, false alarms can be generated with the same magnitude as the true target orientation spectrums. At this point, the FBSS method produces an incorrect orientation output. By applying the method suggested in this paper, it is possible to identify the correct orientations of the targets. However, there is an error between the target orientations obtained and the actual orientations, as can be seen from Fig 6(b). The algorithm's error relies on the spatial spectrum output of the whole array. The advantage is that the spatial spectrum output of the whole array ensures that the algorithm has certain robustness. Nevertheless, it can be seen from the figure that when there is an error in the DOA estimation result, the spatial spectrum output of the array already contains the error. In Fig 6  (c), two pairs of directions in the four targets satisfy the grating lobe position relationship. It can be seen that the proposed method can extract the actual target orientations with low SNR from the grating lobe interferences. Compared with the traditional decoherence FBSS algorithm, the proposed method does not require decoherence operation. In addition, as shown in  The simulation for the existence of the array element position error is carried out to investigate the performance of DOA estimation method in this case. Assume that the array element position error is with a mean of 15% of the spacing between two sensors of the subarray with N sensors. The simulation conditions are the same as Fig 6(c). The effect of SNR on the accuracy of proposed method is shown as Fig 7. Each point is averaged based on 100 trials. One can find that starting from the SNR of -5 dB, the algorithm error gradually stabilizes as SNR becomes larger and the error is about 0.5. As shown in Fig 8, one can find that in the presence of array element errors, the algorithm can still obtain the target directions. Because the algorithm is conducted based on the traditional beamforming and thus inherits the robustness of the traditional beamformer.
In order to further verify the stability of the algorithm, the simulation of the bearing time record (BTR) with low SNR is carried out. At the same time, the algorithm's performance for detecting weak targets under interferences is also investigated. In Fig 9, the parts marked with "SI" represent the strong false alarm interference due to the high grating lobes of the array output. The solid black lines represent the true bearing trails of the sources, and the red asterisks exhibit the DOA estimation results with the proposed method. The interferences sources are from the directions of 48˚and 86˚with SNR = 0 dB, while the true targets come from the directions of 57˚and 78˚with SNR=-5 dB. The target with the azimuth of 48˚is buried in the grating lobe generated by the interference from 78˚and the target with 57˚is buried in the grating lobe generated by the interference of 86˚. Simulation results show that the algorithm can obtain the true orientations of the targets in many interferences. The disadvantage is that a small number of false alarms are inevitably generated due to the strong coherence of the targets.

Conclusion
In this work, we proposed a method for the direction-of-arrival (DOA) estimation of multiple co-frequency targets using a coprime vector sensor array. The problem that many grating lobes with high amplitudes seriously affect the accurate DOA estimation results was investigated. A conventional beamformer for the coprime vector sensor array as a whole was utilized to preprocess the beam data of the array. After obtaining the set of pre-estimated target orientations, the data cross-validation method based on vector cross terms (VCTs) was then used to transform the detection of the co-frequency targets into a solvable processing without spatial smoothing. We presented the Queue Selection (QS) method based on inverse beamforming cancellation technology to reduce the effect of grating lobe interferences further. The accuracy of directions extraction is enhanced by weakening the grating lobes in a sequential interference suppressing. The simulation results show the reliability and stability of the proposed method. Since the algorithm is based on the traditional beamforming output of the array, the disadvantage is that DOA estimation error will be affected by the array output after beamforming. However, it still inherits the robustness of conventional beamformers.